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A numerical study over a nominally two-dimensional circulation control airfoil is performed using a large- 
eddy simulation code and two Reynolds-averaged Navier-Stokes codes. Different Coanda jet blowing condi- 
tions are investigated. In addition to investigating the influence of grid density, a comparison is made between 
incompressible and compressible flow solvers. The incompressible equations are found to yield negligible dif- 
ferences from the compressible equations up to at least a jet exit Mach number of 0.64. The effects of different 
turbulence models are also studied. Models that do not account for streamline curvature effects tend to predict 
jet separation from the Coanda surface too late, and can produce non-physical solutions at high blowing rates. 
Three different turbulence models that account for streamline curvature are compared with each other and 
with large eddy simulation solutions. All three models are found to predict the Coanda jet separation location 
reasonably well, but one of the models predicts specific flow field details near the Coanda surface prior to sep- 
aration much better than the other two. All Reynolds-averaged Navier-Stokes computations produce higher 
circulation than large eddy simulation computations, with different stagnation point location and greater flow 
acceleration around the nose onto the upper surface. The precise reasons for the higher circulation are not 
clear, although it is not solely a function of predicting the jet separation location correctly. 


Nomenclature 


A Planform area 

a Speed of sound 

Cl Section lift coefficient, L/(g 00 A.) 

C p Pressure coefficient, (p — Poo)/qoo 

C rc Constant for SSTRC model 

Cp, Jet momentum coefficient, Eq. (11) 

c Airfoil chord length 

Cbi Constant in SA model 

c r i,c r 2 >c r 3 Constants in SARC model 
c'„ coefficient in EASM-ko model 

Ddi f f Diffusion term in S A model 
Ddiss Dissipation term in SA model 

dist Normal distance from the wall /c 

F\ Function in SST and SSTRC models 

F 4 Sensitization factor, Eq. (9) 

f r i,ft 2 Functions in SA model 

h Slot height 

k Turbulent kinetic energy 

L Lift force 

M Mach number, U/a 

to Mass flow rate 


P Production terms, Eqs. (2) and (8) 

p Pressure 

q Dynamic pressure, pU 1 2 / 2 

Re Reynolds number, pU^c/ p 

Ri Gradient Richardson number, Eq. (10) 

r Radius of curvature 

r Term in SARC model 

r* s/n 

S Magnitude of strain, Eq. (5) 

S Term in SA model 

S,.j Strain rate tensor 

t time 

U Velocity 

u, v, w Cartesian velocity components 

Uj Cartesian velocity component 

ug Streamwise velocity 

uful Specific turbulent stress tensor component 

u'qWq Specific turbulent shear stress, Eq. (12) 

Wij Vorticity tensor 

x, y, z Cartesian coordinates 

Xj Cartesian coordinate 
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Greek: 


Subscript: 

a 

Angle of attack (deg) 

inlet Plenum inlet condition 

(3,(3* 

Constants in SST and SSTRC models 

j Jet condition 

A s+ 

Streamwise spacing in wall units 

max Maximum condition 

An 4 

Normal spacing in wall units 

mean Mean condition 

A6+ 

Spanwise spacing in wall units 

ref Reference (oo) condition 

7 

Constant in SST and SSTRC models 

sep Jet separation location 

P 

Coefficient of viscosity 

oo Free-stream quantity 

Pt 

Turbulent eddy viscosity 


F 

Turbulence variable in SA model 


e 

Angle (deg) around Coanda surface 


p 

Density 


2 

Constants in SST and SSTRC models 


T 

Turbulence time scale 


T ij 

Shear stress 


n 

Magnitude of vorticity, Eq. (4) 


LU 

Specific dissipation rate 



I. Introduction 

Circulation control can dramatically enhance the lift of airfoils and wings. By blowing a jet of air along an airfoil 
surface near its curved trailing edge, separation can be delayed because of the Coanda effect (which causes the jet to 
“stick” to the curved surface). See Fig. 1. Many computational studies have been done for circulation control airfoils; 
for example, see Refs. 1-9. However, in general most efforts have only attained limited success. In particular, many 
Reynolds-averaged Navier-Stokes (RANS) models can do reasonably well for low blowing conditions, but at higher 
blowing rates they tend to predict jet separation too late. In fact, sometimes the computed jet wraps un-physically too 
far onto the airfoil’s lower surface, or even completely around the airfoil. 

There are many difficulties associated with validating codes, methods, and models on circulation control airfoils. 
Due to the typically very narrow blowing-slot openings, the near-wall region of circulation control airfoil experiments 
can be very difficult to measure. In particular, it is important to know the jet exit conditions very accurately, so 
that the CFD can apply the same boundary conditions. The behavior of the jet can be very sensitive to the exit 
conditions. The integral quantity jet momentum coefficient, typically defined from experiment based on total mass 
flow rate and derived jet exit conditions, may be insufficient information for the purpose of CFD validation because 
of possible nonuniform blowing across the entire airfoil span and errors from the isentropic flow assumption. 10 It 
is also particularly difficult to maintain/guarantee two-dimensionality in the wind tunnel flowfield at strong blowing 
conditions. Furthermore, three-dimensional vortical structures in the wind tunnel due to interaction with the side walls 
cause a net downwash (or effective negative angle of attack) that increases in magnitude with increasing blowing. 
Some 3-D RANS studies that included side walls 21011 have demonstrated this effect. The effect of this downwash 
can be difficult to accurately account for when performing 2-D CFD, especially when also attempting to account for 
the upper and lower tunnel walls. To date, 3-D LES computations have typically lacked the resources for attaining 
sufficient grid density to adequately resolve the side-wall vortical structures. 

Of all RANS models applied to circulation control airfoils, arguably the most successful to date have incorporated 
some modeling of curvature effects, which can be significant in circulation control cases near the trailing edge Coanda 
surface. For example, the Spalart-Allmaras for Rotation and Curvature (SARC) model 12 (see Swanson and Rumsey 1 ) 
and a full Reynolds stress model (see Slomski et al. 13 ) tended to yield better results compared to experiment than other 
RANS models like Spalart-Allmaras (SA) 14 or Menter’s shear-stress transport (SST) fc-cu, 15 which do not directly 
account for curvature. However, even curvature-corrected models have not done particularly well at very strong 
blowing conditions; 1 they still tend to yield too high lift compared with experiment. 

Recently, LES computations have also been applied to circulation control flows. There has been some success at 
low blowing conditions, 16-19 but at higher blowing rates LES has not been as successful. 20 The lack of better overall 
success for LES to date may be due to insufficient computational resources. Since computational capacity changes 
very rapidly, this avenue of pursuit is being continued. 

An on-going collaborative experiment at NASA Langley and Georgia Tech has been studying the flow around a 
circulation control airfoil designed at the Georgia Tech Research Institute (GTRI). 21 Because only limited experimental 
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data specifically for the purpose of CFD validation have been obtained to date for this model (see, e.g., Allan et al. 10 ), 
this paper will focus on making comparisons between RANS and LES for flow over the GTRI model, with the intent of 
discerning key differences between the methods. The focus is to employ the methods on identical grids and matching 
boundary conditions at zero angle-of-attack, including the top and bottom walls. Three-dimensional effects due to 
side-walls are neglected, although the LES includes 3-D flowfield structures. We would like to answer the question: 
In what ways do the 2-D RANS and span wise-averaged LES methods differ? Attempting to answer this question 
involves determining the influence of grid density, code, and turbulence modeling. A curvature correction for SST 
(from Hellsten 22 and Mani et al. 23 ), that has not been applied to circulation control airfoils before, is also examined. 
Eventually, it is hoped that LES may be used to improve RANS turbulence models for these kinds of flowfields. 

This paper is organized as follows. After a brief description of the RANS and LES computer codes employed for 
this effort, the grid and test cases are described. Then, in the Results section, four areas of focus are discussed: (a) 
a grid and code sensitivity study, (b) comparisons between compressible and incompressible RANS solutions, (c) the 
effects of turbulence models and curvature corrections for a case with high blowing, and (d) comparisons between 
RANS and LES. Finally, conclusions are drawn. 

II. Numerical Methods and Turbulence Models 

The RANS computer codes employed for this study were FUN3D 24 and CFL3D. 25 FUN3D is an unstructured 
finite volume upwind-biased node-centered code, and CFL3D is a multiblock structured finite volume upwind-biased 
cell-centered code. Both solve the compressible RANS equations and are nominally second order spatially accurate. 

FUN3D solves the equations on mixed element grids, including tetrahedra, pyramids, prisms, and hexahedra and 
also has a two-dimensional path for triangular and quadrilateral grids. It employs an implicit upwind algorithm in 
which the inviscid fluxes are obtained with the flux difference splitting scheme of Roe. 26 For second-order accuracy, 
interface values are obtained by extrapolation of the control volume centroidal values, based on gradients computed 
at the mesh vertices using an unweighted least-squares technique. The solution at each time-step is updated with a 
backward Euler time-differencing scheme. At each time step, the linear system of equations is approximately solved 
with a multi-color point-implicit procedure. Local time-step scaling is employed to accelerate convergence to steady- 
state. For all the applications in this paper, FUN3D solves the turbulence equations with a first-order advection scheme. 

In CFL3D, the convective term is approximated with third-order upwind-biased spatial differencing, and the vis- 
cous terms are discretized with second-order central differencing. The flux difference-splitting method of Roe 26 is 
employed to obtain fluxes at the cell faces. Advancement in time is accomplished via backward Euler, with an implicit 
approximate factorization method. The turbulence equations are solved with a first-order advection scheme. 

The LES code employed is CDR an incompressible Navier-Stokes code developed at the Center for Turbulence 
Research (CTR) at Stanford University. 27 CDP is a finite volume unstructured flow solver with second order spatial 
accuracy. It is based on the least-squares pressure gradient reconstruction of Mahesh et al. 2s for hybrid unstructured 
grids. The methodology is nearly energy conserving in the inviscid limit. CDP uses central differencing and does not 
introduce any numerical dissipation or smoothing. The necessary dissipation comes solely through the subgrid-scale 
model. The subgrid-scale stresses in the CDP are modeled with the dynamic Smagorinsky model. 29 

The RANS turbulence models SA 14 and SST 15 and the LES dynamic Smagorinsky model 29 are standard and will 
not be described here. Details can be found in their respective references. However, a brief discussion will be given 
for the less standard models SARC, EASM-ko, and SSTRC. 

In the SARC model, 12 the production term of the SA model is modified. In the original model, 

dv Ou 

"b 7 = P Ddiff T Ddiss (1) 

where v is related to the kinematic eddy viscosity, P = cvi (1 — fr 2 )Sv, and S' is a function of the magnitude of the 
vorticity. All the terms are defined in the original SA reference. 14 In the SARC modification, 

P = c bl (f rl - f t2 )Sv (2) 

where 

2 r* 

fr l = (1 + c r i)— — - [1 - c I . 3 tan^ 1 (c r2 r)l - c rl (3) 

1 + r* 

with r* = S/Q, and 

|fi| = x 2U „ll 0 . (4) 
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\S\ = x'lS.'Sij. (5) 

and Wij = (du.i/dxj — duj / dxi) /2 is the vorticity tensor and = (dui/dxj + duj/dXi)/2 is the strain-rate tensor. 
The function f depends on the Lagrangian derivative of the strain-rate tensor (see original reference), and the constants 
are taken to be c r \ = 1, c r 2 = 12, and c r 3 = 1. 

The EASM-ko model 30 is an explicit algebraic stress model in k-uj formulation. In this model, the turbulent 
eddy viscosity /j t is given by /j t = c'^pkr, where t = 1 /w for the k-u> version, and c? is a variable coefficient 
obtained through the solution of a cubic equation. This model makes use of a non-linear constitutive relation rather 
than the Boussinesq model. As discussed in Rumsey et ah, 31 because EASM formulations are derived directly from 
the Reynolds stress model, they retain some of the invariance properties of the full differential form even though the 
assumption of anisotropy equilibrium in the Cartesian frame of reference is not strictly correct for curved flows. A 
curvature corrected form of EASM was derived, but curved flow results were very similar to those predicted by the 
standard form of the model. Because the standard model is more robust than the curvature corrected version, only 
EASM-ko (the standard model) was used in the current work. 

SSTRC is Hellsten’s curvature correction 22 to the standard SST model. In SSTRC, the destruction term in the 
w-equation is multiplied by a curvature sensitization factor F 4 , which is a function of the local flow curvature. The 
equations are: 


dk dk d 

p m +plll ar, = ^ p aT, , 


M- 


Pt\ dk 


a k J dx 


3 J 


du , du 7P D P o 2 . d 

•‘si + a ^ + &; 


Pt \ du 


dxj 


2 P 


1 — F 1 dk duj 
dx-j dxn 


where the turbulence production term is 


P ~ Ti °dx 


du^ 
3 


The sensitization factor is given by 


(6) 

(7) 

( 8 ) 


Fa 


1 

1 + G rc Ri 


(9) 


The curvature sensitization factor F 4 is a local formulation that is independent of the coordinate system and based on 
Hellsten’s approximation to the gradient Richardson number: 


Rl \s\\\s\ 7 ’ (10) 

This gradient Richardson number is easily computed, and is independent of coordinate system and grid topology. 
Near-wall convex curvature yields F 4 levels below 1, which decreases ^.’-destruction and in turn increases to. The 
near-wall eddy viscosity decreases as a result. Therefore, usage of the correction coefficient F,\ weakens the turbulent 
boundary layer around convex curved surfaces, which creates a greater tendency to separate. 

Hellsten originally recommended G rc = 3.6, but this value was found to diffuse the jet excessively in some 
applications. Mani et al. 23 recommended a value of C rc = 1.4 for cases involving ground-jet interactions. This value 
is also used for all results here. See the original references for details concerning other terms in the equations. 


III. Description of Grid and Test Cases 

A near view of the 2-D grid (with every other gridpoint removed) is shown in Fig. 2. This figure also shows the 
geometry of the airfoil, including the dual-chambered plenum. The earlier sketch of Fig. 1 showed the geometry near 
the Coanda surface (note that the lower surface included a slot for blowing, but it was closed for the current study). 
Although not shown, the top and bottom tunnel walls were located at z = ±0.508 m (z/c = ±2.327). A wavy 
transition strip was located near the lower surface leading edge (for details, see Nishino et al. 17 ). The airfoil chord 
was c = 0.2183 m. The jet slot height was h = 0.000503 m ( h/c = 0.0023). The radius of the Coanda surface was 
r = 0.02066 m (r/c = 0.09463). 

The fine grid employed for the LES computations had approximately 116 million grid points. It had a spanwise 
(y-direction) extent of 0.014 m (or y/c = 0.0641), with 256 points in that direction and periodic boundary conditions 
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at the side boundaries. For the lowest jet blowing condition, the coarsest grid resolution (in wall units) over the 
Coanda surface was approximately As + « 42, A n + « 0.6, A/; ss 25, in the streamwise, wall-normal, and spanwise 
directions, respectively. Note that the resolution in wall units reported in Table 1 of Nishino et al. 17 was not correct. 
The subgrid-scale eddy viscosity in the attached jet region (i.e., before the separation of the jet from the Coanda 
surface) was up to about 5 times larger than the molecular viscosity. For the higher jet blowing condition, the coarsest 
resolution was approximately As + « 60, A n + « 0.8, A b + « 36, and the subgrid-scale eddy viscosity in the attached 
jet region was up to about 10 times larger than the molecular viscosity. The grid included the tunnel top and bottom 
walls, and also modeled the two inner plenum chambers as well as the wavy transition strip on the lower surface near 
the leading edge. The fine grid used for the RANS computations was the same, but was solved two-dimensionally 
using 451,400 grid cells in one spanwise row. It is important to note that for the 2-D RANS computations, the wavy 
strip (described in detail in Nishino et al. 17 ), is simply a small rectangular bump located near x/c= 0.021 with height 
0.00256c and length 0.01368c. The grid was originally constructed as a multi-zone grid with 73 zones. For use with 
FUN3D, a single zone consisting of hexahedral elements was created from the original grid by removing all duplicate 
points at the zone interfaces. 

The free-stream flow conditions were: M = 0.099, Re = 2.24 x 10 6 per m, or Re = 488,992 per chord. The 
RANS computations were run in steady-state mode, whereas the LES computations were run time-accurately and 
averaged to obtain statistically steady solutions, as described in Nishino et al. 17 Because the simulations included 
wind tunnel walls, no angle of attack corrections were applied to any cases; all were run at a = 0°. Viscous wall 
boundary conditions were applied on all surfaces of the airfoil, including the interior plenum. The only exception was 
the left-most face of the interior-most plenum, at which a velocity and density were specified (this face had a height 
of 0.01261 m). The top and bottom tunnel walls were treated as inviscid surfaces. The upstream inflow plane (located 
at x/c = —3.621) had a farfield Riemann invariant boundary condition, and the downstream outflow plane (located at 
x/c = 10.005) had a specified back pressure p/p re f = 1-0. 

The various plenum boundary conditions employed are listed in Table 1 . Conditions were chosen to match the 
mass flow rate, with total temperature at inflow of 295 K. The jet momentum coefficient is defined as: 


C,= 




mean 


QrefA 


(ID 


where rhj is the mass flow rate at the jet exit and Uj tmean is the mean jet velocity at the jet exit. Note that the specific 
conditions for compressible and incompressible boundary conditions at the inlet for the RANS and the LES needed to 
be somewhat different in order to achieve a similar velocity at the jet exit. The incompressible RANS computations 
in this paper used the (A) conditions; however, the difference between results using (A) and (B) conditions was fairly 
small. The reference velocity was U re f = 34 m/s. The reference speed of sound was a re f = 344 m/s and the 
reference density was p re f = 1.2 kg/m 3 . Most comparisons were performed for the two lowest blowing conditions 
( Mj « 0.39 and 0.64). Only a few compressible RANS computations (and no LES) were performed at the highest 
blowing rate of Mj ss 0.90, in order to demonstrate the effect that the various turbulence models have at that high 
blowing rate. 

It is recognized that the local Mach number of the flow near the jet exit can exceed levels for which incompressible 
flow assumptions are considered valid. However, at least for the two lower blowing cases, it was felt that the high- 
Mach-number regions were limited to a fairly small area so that use of the incompressible LES code was justified. 17 As 
will be shown below, FUN3D was run both in compressible and incompressible modes and its results were compared 
in order to quantify the influence of solving the incompressible equations and justify their use for these flows. 


IV. Results 

A list of the CFD computations performed is provided in Table 2. The table includes computed lift coefficients 
and jet separation angle ( 9 sep ), as defined in Fig. 1. Note that the Mj = 0.39 jet separation angle for the LES listed 
here is slightly different from that given in Nishino et al. 18 because here we are accounting for the very thin region of 
mean reverse flow near the body that occurs near the jet separation point. In other words, the separation angle for all 
results in this paper is defined as the location where wall skin friction goes to zero in the time-averaged flow field. In 
the “extra-fine” grid used for one of the LES computations, the number of tangential cells around the Coanda surface 
was doubled (from 800 to 1600). All other aspects of the grid remained the same as the fine grid. On the other hand, 
the “medium” grid level used for the CFL3D grid-independence studies consisted of every other point removed from 
the fine grid in each coordinate direction. 

Fig. 3 shows time-averaged surface pressure coefficients from the LES computations, demonstrating that - other 
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Table 1. Summary of Coanda blowing conditions 


conditions per unit span 




Uinlet-, m/ S 

Piniet, kg/m 3 

to = 0.0636kg/(m s) (C M s= 

a 0.044, Mj ? 

« 0.39, Uj^max ? 

a 130 m/s) 



compressible RANS 




3.925 

1.285 

incompressible RANS (A) 




4.268 

1.182 

incompressible RANS (B) 




4.203 

1.2 

incompressible LES 




4.203 

1.2 

to = 0.1048kg/ (m s) (C p s 

a 0.120, Mj ? 

^ 0.64, U j^rnax ^ 

a 210 m/s) 



compressible RANS 




5.703 

1.457 

incompressible RANS (A) 




7.126 

1.166 

incompressible RANS (B) 




6.924 

1.2 

incompressible LES 




6.924 

1.2 

to = 0.1403kg/(m s) (C M s 
compressible RANS 

a 0.215, Mj ? 

- 0 QQ TJ £ 

^ kj.ukj, kj j^max n 

a 300 m/s) 

6.248 

1.781 


than in the immediate region of the wavy transition strip - the results were essentially statistically two-dimensional (see 
also Nishino et al. 17 ). Therefore, direct comparisons between spanwise-averaged LES and 2-D RANS were possible. 

A. Grid and Code Sensitivity Study 

LES computations were performed for the Mj = 0.64 case using both the fine grid (approximately 116 million total 
grid points) and the extra fine grid with double the number of tangential cells around the Coanda surface (approxi- 
mately 175 million total grid points). Results were very nearly the same. Visualizations of the flow structures around 
the Coanda surface as well as averaged surface pressure coefficients are shown in Fig. 4. As described in Nishino et 
al., 17 a number of hairpin-like flow structures were observed in the LES around the mixing layer of the high-speed jet 
flow and the low-speed external flow. Here, these flow structures appeared to be similar between the two grids. 

Many 2-D RANS computations were performed on both the fine grid (45 1 ,400 grid cells) as well as on the medium 
grid with every other point removed ( 1 12,850 grid cells). Furthermore, the two codes FUN3D and CFL3D were run on 
the same fine grid using the common SA, SST, and SSTRC turbulence models as a code-independence check. FUN3D 
did not have the ability to run the SARC or EASM-ko turbulence models, but these were run with CFL3D on the 
two grid sizes. Results are tabulated in Table 2, and surface pressure coefficients are shown in Fig. 5. For the most 
part, results were nearly the same on the two grid sizes and with the two codes. (However, note that EASM-ko did not 
produce a steady-state solution at Mj = 0.39 on the medium grid, so only the fine grid result is given.) Two exceptions 
to this were the SA and SST models on the Mj = 0.90 case. These particular solutions were non-physical and not 
well-converged: the SA model allowed the jet to completely wrap around the airfoil, and the SST model allowed the 
jet to stick to the Coanda surface all the way to the lower surface, so that the jet faced forward into the oncoming 
stream. In the c p plots, this non-physical behavior shows up as a large “kink” or “step” in the lower surface pressure 
coefficient between x/c = 0.90 and 1.0. This behavior was consistent between codes and occurred on both fine and 
medium grids. 

Overall, the differences between the turbulence models were much greater than the differences due to grid or code. 
This fact is important, because it makes it possible to draw conclusions about the behavior of the turbulence models 
with some degree of confidence that discretization errors and code bias are probably not playing a significant part. As 
seen in Table 2, the effect of grid/code on lift coefficient and jet separation location were also minimal. 

B. Comparisons Between Compressible and Incompressible RANS Solutions 

For even moderate blowing conditions, the jet behaves compressibly near the jet exit. Nonetheless, researchers have 
employed both incompressible and compressible CFD codes to circulation control airfoil flows, but until now the 
impact of the type of solver has not been explored. It can be argued that the compressible region is very localized for 
many blowing rates, and hence the use of incompressible equations may not have much negative consequence for the 
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Table 2. Summary of computations performed 


Method 

Mj 

Code 

Grids 

Turbulence model 

c L 

Jet separ. ( 0 sep ) 

Incompressible LES 

0.39 

CDP 

fine 

dynamic Smag 

1.36 

67° 


0.64 

CDP 

fine / extra-fine 

dynamic Smag 

3.50/3.51 

95° / 96° 

Incompressible RANS 

0.39 

FUN3D 

fine 

SA 

1.85 

73° 


0.64 

FUN3D 

fine 

SA 

4.67 

107° 

Compressible RANS 

0.39 

FUN3D 

fine 

SA 

1.84 

72° 


0.39 

FUN3D 

fine 

SST 

1.67 

68° 


0.39 

FUN3D 

fine 

SSTRC 

1.56 

64° 


0.64 

FUN3D 

fine 

SA 

4.65 

107° 


0.64 

FUN3D 

fine 

SST 

4.49 

106° 


0.64 

FUN3D 

fine 

SSTRC 

3.93 

95° 


0.90 

FUN3D 

fine 

SA 

4.91 

wrap 


0.90 

FUN3D 

fine 

SST 

6.46 

174° 


0.90 

FUN3D 

fine 

SSTRC 

6.69 

117° 


0.39 

CFL3D 

fine / med 

SA 

1.86/1.82 

73° / 73° 


0.39 

CFL3D 

fine / med 

SARC 

1.72/1.71 

68° / 67° 


0.39 

CFL3D 

fine / med 

SST 

1.70/1.69 

69° / 68° 


0.39 

CFL3D 

fine / med 

SSTRC 

1.59/1.57 

65° / 64° 


0.39 

CFL3D 

fine / med 

EASM-ko 

1.70/unst 

68° / unst 


0.64 

CFL3D 

fine / med 

SA 

4.72/4.66 

108° / 108° 


0.64 

CFL3D 

fine / med 

SARC 

4.04 / 4.02 

95° / 95° 


0.64 

CFL3D 

fine / med 

SST 

4.42 / 4.40 

105° / 103° 


0.64 

CFL3D 

fine / med 

SSTRC 

3.92/3.90 

95° / 93° 


0.64 

CFL3D 

fine / med 

EASM-ko 

4.10/4.07 

97° / 96° 


0.90 

CFL3D 

fine / med 

SA 

4.61/4.50 

wrap / wrap 


0.90 

CFL3D 

fine / med 

SARC 

6.96/6.84 

112° / 112° 


0.90 

CFL3D 

fine / med 

SST 

6.88/7.02 

169° / 155° 


0.90 

CFL3D 

fine / med 

SSTRC 

6.79/6.69 

119° / 118° 


0.90 

CFL3D 

fine / med 

EASM-ko 

6.74/6.71 

116°/ 115° 
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global flow field. 

By using a single code (FUN3D) that has both incompressible and compressible capabilities, the impact is explored 
for the first time here. Figs. 6(a) and (b) show streamlines for the Mj = 0.64 case, using FUN3D in incompressible and 
compressible mode, respectively. Fig. 6(c) shows the streamwise velocity (ug nondimensionalized by the reference 
velocity) profile comparisons. Incompressible and compressible results were essentially the same, including the jet 
separation location. Although not shown, results for the Mj = 0.39 case were even closer, as would be expected. The 
Mj = 0.90 case produced non-physical results using the SA model, so did not provide any useful information. These 
results appear to justify the use of incompressible flow solvers (like the CDP code employed for LES computations), 
at least for blowing rates up to Mj = 0.64. 

C. Effect of Turbulence Model and Curvature Correction for Mj = 0.90 

Results from Swanson et al. 2 indicated that a rotation/curvature correction to the SA turbulence model can signifi- 
cantly improve results for circulation control cases. Furthermore, Swanson et al. showed that the SST model tended 
to produce poor solutions at higher blowing rates. Here, two other turbulence models that provide some level of sen- 
sitivity to flow curvature - EASM-ko and SSTRC - are also examined for the Mj = 0.90 case. EASM models were 
also investigated previously for circulation control airfoils by Gross and Fasel. 8 9 

Earlier, Fig. 5 showed surface pressure coefficients for all the RANS models on all three cases. The two cases that 
stood out (in a negative sense) were SA and SST at Mj = 0.90. Figs. 7(a)-(e) show streamlines for the Mj = 0.90 
case, using all five RANS models. The jet in the SA result wrapped completely around the airfoil (Fig. 7(a)), and 
the jet in the SST result separated close to the end of the Coanda surface near 6 = 170° (Fig. 7(b)). Neither of 
these models accounts for curvature. The three models that account for curvature yielded separation locations near 
6 « 110° — 120°. Based on other circulation control airfoil experiments (e.g., Novak et al. 32 ), this separation location 
appears much more reasonable than that predicted by SA or SST. 

D. Comparisons Between RANS and LES Solutions 

The mean LES lift coefficients as a function of C ; , are compared with the steady 2-D RANS results in Fig. 8. Generally, 
the SARC, SSTRC, and EASM-ko results were within a few percent of each other at all blowing values. The SA and 
SST levels were higher at the Mj = 0.64 condition. At the highest blowing condition, the SA wrap-around produced 
a very low Cl value, while the SST (coincidentally) produced Cl in agreement with SARC, SSTRC, and EASM-ko 
in spite of the fact that its jet separation angle was so much larger. The LES results were lower than the best of the 
other models by between 12 and 17%. 


1. Mj = 0.39 Case 

Fig. 9 shows surface pressure coefficients for the Mj = 0.39 case. As discussed in Nishino et al., 17 the LES yielded 
excellent agreement with (preliminary) experimental data at this blowing condition. The LES also agreed well with 
independent implicit LES computations of Madavan and Rogers 19 (note that the latter computations were not direct 
numerical simulations, as implied by the paper’s title). Here, all RANS results predicted more upper surface suction 
than the LES. The jet in the SA model result separated slightly later (by 5-6°) than LES, but all the other RANS models 
predicted the separation within 1-3° of LES. 

Both here and in Nishino et al., 18 C. p results using SA were furthest from LES. Here, the SST, SARC, and EASM- 
ko were all about the same, and SSTRC was the closest to LES. Additional details are shown in Figs. 10 and 11. (In 
these figures, only the best RANS results SSTRC, SARC, and EASM-ko are shown.) The SSTRC and EASM-ko 
models both predicted a higher peak jet velocity and narrower jet width compared to LES, while the SARC model 
was much closer. In fact, SARC and LES velocity profiles were very similar prior to jet separation. In Fig. 11, the 
nondimensional turbulent quantity u' g w' e is given in a coordinate system aligned with the local near-wall flow and 
flow-normal directions, where 


u e w e = 



u'u')sin{26 ) + u'w'cos{26 ) 


( 12 ) 


with u'u ' , w'w', and u'w' the Cartesian specific turbulent shear stress components. Note that the LES results in- 
cluded only the resolved component (i.e., the subgrid-scale contributions, which were fairly small in the attached jet 
region, were not included). On the whole, all RANS models overpredicted the u' s w' g peak magnitude close to the wall 
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compared to LES. In the outer part of the jet, the SSTRC and EASM-ko models severely underpredicted the peak 
magnitude and extent, whereas the SARC model was closer in the regions prior to jet separation. 

Comparing the streamlines around the back of the airfoil in Fig. 12, it can be seen that SSTRC yielded a bubble 
that was much too large compared to LES, and EASM-ko yielded a separation bubble that was very different in 
character from the LES. The SARC solution was very close to the LES solution, but the separation bubbles were 
slightly different in shape, and the SARC streamlines drove vertically downward further than for LES. Fig. 13 shows 
jet/wake profiles behind the airfoil. At the first position near the closure of the separation bubble, all the RANS models 
were fairly consistent with LES, but further downstream there were bigger differences. The RANS models all almost 
completely dissipated the very weak jet, and exhibited more wake -like behavior than LES. With these differences in 
bubble size/position and jet/wake behavior and position, the RANS solutions had a resulting larger overall circulation 
than LES. Thus, as shown in Fig. 14, the RANS stagnation point near the airfoil leading edge was further downstream 
(note that the region near the trip is not shown in the LES solution because the 3-D spanwise waviness of the trip 
itself precluded spanwise averaging in this area). Although not shown, the other RANS models were similar to SARC. 
Therefore the RANS flow accelerated more around the nose. Boundary layer profiles on the upper and lower surfaces 
near x/c = 0.8 are shown in Fig. 15. Both upper and lower surfaces were turbulent at these locations. Because the 
stagnation point locations did not agree, the edge velocities for RANS were different than for LES, as expected. 

2. Mj = 0.64 Case 

The Mj = 0.64 surface pressure coefficients and Coanda profile results are shown in Figs. 16, 17, and 18. Similar 
to the lower blowing case, here the SA and SST yielded the worst agreement with LES. The SA and SST produced 
similar jet separation locations between 105-108°, whereas LES, SSTRC, SARC, and EASM-ko all separated earlier, 
near 95-97°. Among the RANS models that accounted for curvature, all produced very similar c p distributions, but 
none agreed with LES. Again, SARC came the closest to agreeing with the velocity and u' g w' g profiles from LES prior 
to separation. 

The Mj = 0.64 blowing case was different in character than the Mj = 0.39 case, because there was no longer a 
large separation bubble directly next to the underside of the separating jet, but rather a much smaller one further around, 
nearer to the bottom of the Coanda surface. (As was seen in Fig. 7, at the even stronger blowing with Mj = 0.90 
this separation bubble went away completely.) Fig. 19 shows details of the flowfields near the trailing edge. The 
separation bubble was slightly smaller for SARC than for LES, but SSTRC and EASM-ko agreed well with LES in 
this regard. In spite of the fact that the RANS solutions separated at about the same location as the LES solution, 
the RANS jets moved downward further vertically than LES downstream. As seen in Fig. 20, in this higher blowing 
case the jet remained fairly strong post-separation, and in terms of peak velocity and spreading rate, the SSTRC and 
EASM-ko both did a reasonably good job predicting it. However, the SARC model dissipated the jet too quickly. This 
is a known problem with the SA model itself, 33 and the problem no doubt carries over to the SARC variant. Again, 
all RANS models produced more circulation than LES. As seen in Fig. 21, the leading edge stagnation point for the 
RANS model was again located further aft on the lower surface, causing higher velocities as the flow accelerated onto 
the upper surface (although not shown, the leading edge behavior of SSTRC and EASM-ko were similar to SARC). 
Boundary layer profiles on the upper and lower surfaces near x/c = 0.8 in Fig. 22 indicate that both upper and 
lower surfaces were turbulent at these locations. Again, because the stagnation point locations did not agree, the edge 
velocities for RANS were different than for LES. 

V. Discussion and Conclusions 

In past studies, it has been generally assumed that if one could capture the correct jet separation location from the 
Coanda surface of circulation control airfoils, then the simulation would yield accurate airfoil forces and pressure dis- 
tributions. Many “standard” turbulence models in earlier studies tended to predict jet separation too late, consequently 
yielding lift predictions that were too high (compared to experiment) by a significant amount. Here, although detailed 
experimental data (including measurements of jet boundary conditions) for the particular configuration being studied 
have yet to be completed, we have compared 2-D RANS computations using various turbulence models with 3-D LES. 
We have demonstrated that turbulence models that account for curvature can produce jet separation locations that agree 
very well with LES. Three different models - SSTRC, SARC, and EASM-ko - all performed about the same in terms 
of the jet separation location and surface pressure coefficient predictions. However, the SARC model yielded much 
better detailed agreement in terms of velocity and turbulence profiles over the Coanda surface compared to LES. Past 
studies (e.g., Swanson and Rumsey, 1 Pfingsten et al. 7 ) have also shown the SARC model to be reasonably good for 
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these types of flows. 

Yet, in spite of agreeing fairly well for many details of the flowfield over the Coanda surface, the best RANS 
models still produced larger airfoil circulation (higher lift) than LES, by between about 12 and 17%. The specific 
reasons for the larger circulation have not been determined yet. Some possibilities include: (a) differences in the shape 
of the separation bubble underneath the jet separation and the direction and location of the separated jet sheet, (b) 
inadequate spreading or dissipation of the downstream jet in the RANS models, and (c) insufficient grid resolution for 
LES, RANS (or both) in regions farther from the airfoil. Also, although SARC was the best at predicting the flow 
very near the Coanda surface, it was worse than the SSTRC and EASM-ko models in predicting the free-shear jet 
characteristics downstream in the M 3 = 0.64 case. The turbulence model must be able to accurately capture both 
near- wall curvature effects and free-shear planar jet decay for these problems. It appears that without a sharp trailing 
edge (Kutta condition), the jet behavior downstream can be extremely sensitive in these circulation control airfoil 
cases. And the airfoil circulation in turn appears to be very sensitive to small differences in jet behavior and location. 

In terms of grid resolution, the RANS made use of one plane of the 3-D LES grid; this was very fine compared 
to grid sizes typically used for 2-D RANS in past studies. It was demonstrated that this grid size was fine enough 
to yield RANS solutions with very little differences from solutions on half the grid size in each coordinate direction. 
Two independent RANS codes also produced almost the same results. It is not known at this time whether the LES 
computations were adequately grid-converged. Best-practices for airfoil-type LES computations were followed within 
the available resource constraints. One LES computation required on the order of 2 months to complete using 256 
processors. (Each fine grid RANS computation took on the order of several hours on 6 processors). A refinement 
of the grid in the region of the Coanda surface made little difference to the LES solution, but only the streamwise 
direction was refined due to resource constraints. 

In all of the computations performed, care was taken to set the boundary conditions inside the jet plenum in order 
to match both the jet mass flow rate and the velocity profile at the jet exit. With this same starting condition at the 
jet exit for all codes, it was possible to confidently conclude that differences downstream were due to the turbulence 
modeling, and not due to other factors such as boundary condition influence. 

As mentioned above, detailed experimental data for this configuration for the purpose of CFD validation is being 
acquired, but is not yet available. Preliminary comparisons between LES and experiment 17 indicated near perfect 
agreement between LES and experiment at the lowest blowing condition, but poorer agreement at a higher blowing 
condition. Measured experimental jet boundary conditions at the slot exit will be required before drawing any firm 
conclusions, however. Also, it is well-known that interactions of the jet sheet with tunnel side walls cause roll-up 
vortices downstream, which induce a net downwash on the airfoil. This effect increases with increasing blowing rate. 
To date this downwash influence has not been accounted for in the LES computations. 

In conclusion, this paper has highlighted some key differences between RANS and LES for circulation control 
airfoil flow. It was shown that use of an incompressible LES code was likely to be valid for this configuration for jet 
blowing up to at least Mj = 0.64. This was demonstrated by using an incompressible and compressible version of the 
same RANS code, and showing essentially identical results. In other words, for these cases the small region near the jet 
exit where the local Mach number exceeded levels where the incompressible assumption is strictly valid did not affect 
the solution noticeably. Several RANS turbulence models were compared with each other and with LES. At very high 
blowing rates, the standard SA and SST models, which do not account for flow curvature effects, had a tendency to 
produce non-physical solutions, with SA somewhat worse than SST. Three different models that account for curvature 
- SSTRC, SARC, and EASM-ko - all behaved better, and were generally similar to each other overall. The behavior 
of SSTRC on a circulation control airfoil case has not been documented previously. All three models did a good job 
predicting the jet separation location. SARC agreed much better with LES in terms of Coanda profile details, but 
worse in terms of downstream jet dissipation and spreading. All RANS models produced higher circulation (higher 
lift) than LES, with different stagnation point location and greater flow acceleration around the nose onto the upper 
surface. The precise reasons for the higher circulation are not clear, although it is not solely a function of predicting 
the jet separation location correctly, as previously believed. A current hypothesis is that differences in the predicted 
character and location of the jet sheet after it separates from the airfoil are a significant contributing factor. 
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Figure 1. Sketch of Coanda flow. 
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Figure 2. View of circulation control airfoil grid (with every other grid point removed). 
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Figure 3. Time-averaged surface pressure coefficients from LES computations. 
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(a) Iso-surfaces of second invariant of velocity gradient ten- 
sor, colored by the magnitude of velocity (uiUi ) 1 ' 2 /U re f, 
fine grid 



(b) Iso-surfaces of second invariant of velocity gradient ten- 
sor, colored by the magnitude of velocity (uiUi ) 1 ' 2 /U re f, 
extra-fine grid 



(c) Surface pressure coefficients 


Figure 4. Effect of grid refinement on LES results, Mj = 0.64. 
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Figure 5. Effect of code and grid refinement on RANS surface pressure coefficients. 
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(a) Incompressible streamlines superimposed on eddy vis- 
cosity contours (arrow indicates jet separation location) 
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(b) Compressible streamlines superimposed on eddy vis- 
cosity contours (arrow indicates jet separation location) 
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(c) Velocity profiles 

Figure 6. Comparisons of incompressible and compressible results for Mj = 0.64, FUN3D, SA turbulence model. 
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Figure 7. Streamlines (superimposed on Mach contours) for compressible RANS simulations, Mj = 0.90, CFL3D. 
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Figure 8. Lift coefficients as a function of blowing coefficient. 



Figure 9. Surface pressure coefficients, Mj = 0.39. 
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Figure 10. Velocity profiles, Mj = 0.39. 
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Figure 11. Turbulent u' e w' e profiles, Mj = 0.39. 
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Figure 12. Streamlines near trailing edge (superimposed on velocity magnitude contours), Mj = 0.39. 
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Figure 13. Comparison of jet/wake profiles, Mj = 0.39. 
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Figure 14. Comparison of LES and RANS velocity magnitude contours near leading edge, Mj = 0.39. 
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(a) Upper surface, velocity profiles 


(b) Upper surface, turbulent profiles 
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(d) Lower surface, turbulent profiles 


Figure 15. Boundary layer profiles near x/c = 0.8 on upper and lower surfaces, Mj = 0.39. 
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Figure 16. Surface pressure coefficients, Mj = 0.64. 
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Figure 17. Velocity profiles, Mj = 0.64. 
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Figure 18. Turbulent u' e w' e profiles, Mj = 0.64. 
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Figure 19. Streamlines near trailing edge (superimposed on velocity magnitude contours), Mj = 0.64. 
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Figure 20. Comparison of jet/wake profiles, Mj = 0.64. 
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Figure 21. Comparison of LES and RANS velocity magnitude contours near leading edge, Mj = 0.64. 
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Figure 22. Boundary layer profiles near x/c = 0.8 on upper and lower surfaces, Mj = 0.64. 
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